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Abstract 



We prove strong convergence of conforming finite element approximations to the stationary Joule 
heating problem with mixed boundary conditions on Lipschitz domains in three spatial dimensions. 
We shov^f optimal global regularity estimates on creased domains and prove a priori and a posteriori 
•^ ' bounds for shape regular meshes. 
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)Q ■ 1 Introduction 

O ■ 
en . 

\^ , The stationary Joule heating problem is a two way coupled system of non-linear elliptic partial differen- 

-^ ■ tial equations modelling the heat and electrical potential in a body. The electrical current acts as a heat 

O ' source in a resistive material while the temperature feeds back to the electrical potential through the 

^ i electrical conductivity. Joule heating is important in many micro-electromechanical systems, where the 

effect is used to achieve very exact positioning at the micro scale, e.g. [15] . The Joule heating problem is 
also studied for the design of semiconductors, in particular in the setting of thermistors. In applications 
/\ . boundary conditions of mixed type are typically used. 

S I The main difficulty in proving the existence of finite energy solutions to the Joule heating problem is 

that, given a finite energy potential, the source term of the heat equation is in general only in L^ , which 
means that the usual variational framework is not directly available. This issue has been studied in [7, 13, 
4] , for Dirichlet boundary conditions, and later in [8, 17] , for mixed boundary conditions. IVIultiplicity of 
solutions and stability was studied in [9] . Similar questions have also been raised for the time dependent 
case, see for instance [2, 27, 22]. 

There have been several works on the numerical solution of the Joule heating problems in recent 
years. For the steady state formulation both conforming and non-conforming finite element methods 
have been studied using homogeneous Dirichlet boundary conditions [30, 31] and homogeneous mixed 
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boundary conditions [29] . Under assumption of sufficient regularity of the solution and small data, a pri- 
ori error bounds have been derived with convergence rates. There have also been parallel investigations 
into numerical methods for the time dependent Joule heating problem, see e.g. [11, 1]. The assumption 
on small data can here be avoided since the Gronwall lemma is available. An a posteriori error: bound for 
a time dependent obstacle thermistor problem is presented in [3] . 

In this paper we prove the strong convergence (of subsequences in case of non-unique exact solu- 
tions) of Galerkin approximations to finite energy solutions of the Joule heating problem in three dimen- 
sions with mixed boundary conditions, using only very mild assumptions on the computational domain 
and the data. The analysis covers, in particular, conforming finite element approximations. To this end 
we introduce a truncation operator in the approximate potential without affecting the consistency of the 
method. Thereby we avoid the assumption of L°° bounds on the discrete potential solution, indepen- 
dent of the mesh size, which are used in [16]. These L°° bounds are very difficult to realize in practice 
in three spatial dimensions. They also impose restrictions on the computational meshes as well as the 
order of convergence of the method. Under the assumption of a so-called creased domain together with 
a sufficiently weak temperature dependency in the electrical conductivity we also prove optimal global 
regularity estimates together with local estimates guaranteeing smooth solutions away from the bound- 
ary given smooth data. We further prove a priori and a posteriori error bounds for conforming finite 
element approximations on shape regular meshes. In our analysis the small data assumption relaxes as 
the coupling of the equations weakens. 

The paper is organized as follows. In Section 2 we introduce the stationary Joule heating problem 
with mixed boundary conditions. In Section 3 we study the continuity properties of the differential op- 
erator to show the convergence of Galerkin approximations to finite energy solutions. In Section 4 study 
the global and interior regularity of solutions on creased domains. Finally, in Section 5 we derive optimal 
a priori and a posteriori error bounds for h-adaptive conforming finite element approximations to the 
Joule heating problem. 

2 The Statement of the Stationary Problem 

Let f2 be a bounded Lipschitz domain in IR^. Let D^ and Du be subsets of 3Q, whose boundaries SD^ and 
dDu are Lipschitz regular manifolds and set iV^ := dD. \ D^p and Ru := dCl \ Du- We shall impose Dirichlet 
boundary conditions for cp and u on D^ and Du, Neumann conditions for cp on A/,^ and Robin conditions 
for u on Ru- 

The scale of Sobolev spaces is denoted by wf. Let, for 5 > 1/p, 

WP(Q;D^) := {v e WfCQ) : v\d^ = 0}. 

Define W^[D.;Du) analogously and denote W^ spaces also with H^. 

Assume for the Dirichlet data that g^ e Wf (Q) n L°° (H) and gu e W^^^ (^^) ^^^ for the Robin data that 
h e W^^j^^dn,). Let a e C^ [U) be bounded from below by a positive Co e K and from above by (t° eU and 
let K e L°°(_Ru) be non-negative. Assume that there are the Poincare-Friedrichs inequalities 



\LHa)< \\yVhHa.m Vi//e VKf (Q;D^), 

ll"^llL2(n)^ IIVw||i2(s^^K3) + ||v/Kif||i2(^j VweWf{n;Du). 

Allow Du = provided the Poincare-Friedrichs inequality remains valid. 



(1) 



The strong formulation of the Joule heating problem is to find (f) e W^(Q) and u e W^{D.) such that 

-V-{o-(w)V(/))=0 o -a[u)Au = ya[u)-V(f), (2) 

-A w-cr(w) I V(/)|2=0 

with the Dirichlet conditions cplu^ = g^p, u\u„ = gu and the natural boundary conditions dyCp = on N(f, 
and KU + dyU = h on Ru with the outward unit normal v. 

Remark 1. For some applications the Lipschitz assumption on Q is too restrictive, a good examples being 
geometries which locally resemble the two-brick domain. We point out that Theorem 1 remains valid for 
domains for which the usual Sobolev embedding holds, a trace operator is available and integration-by- 
parts can be carried out. For example, see [20] for more information in this direction. For Theorem 2, 
however, the Lipschitz assumption is an essential part of the definition of creased domains. 

2. 1 The Weak Formulation of the Stationary Problem 

A weak solution of the stationary Joule heating problem is a 

(0, u) = ig^ + (p, gu + meH^ m X //I (o) 

such that £ H^{Q.,D^), u e H^{Q.,Du) and 

<o-(w)V0,Vi//> = 0, ] 

\ (3) 

{Vu.Vw) + {ku,w)r^^ = {a[u)'\I(p-y(p,w) + {h,w)R^^ J 

foralli//e H^[p.;D(p) and V we W^{p.;Du). Indeed, the choice of spaces ensures that cjCmjV^-V^ eL^m) 
which guarantees that the second equation is meaningful for all w e W^iP.; Du). 

Lemma 1. If{(p, u) is a solution of (3) then 

go < (/) < g°, with g° = maxg0, go = min g^. 

xeD^ x€.D^ 

Proof Define J = max(0,(/)-g°) e H^{VL;D^). One can use j as a test function in equation (3): 

= <o-{M)V(/), Vx> = <o-{M)V(</.-g°), Vx> 

= /suppwnn «■(") ^I • Vj dx = <o-(w) Vj, Vj) . 

Now use the Poincare-Friedrichs inequality to get llxllL2(n) = 0, so (/) < g°. An analogous argument with 
gogives(/)>go. D 

Because of the maximum principle we may introduce an equivalent weak formulation which em- 
ploys the cut-off functional 

[/I := min(max(/ + g^, go), g°) - g^. 

Then go - g^ < [/I < g° - g^f, and \(p] = (p. This functional is essential in the proof of the convergence 
of Galerkin solutions without the need for a discrete maximum principle; a property desirable from the 
numerical point of view. 



Lemma 2. The set of functions which satisfy 

<o-{m)V(/),Vi//> = 0, 

{yu,yu>} + {Ku,u>}R,, = > (4) 

-{a[u) \(f)]V(f),Vw) + {a{u)Vg^-V(l), w) + {h, w)r^, 
for all (i//, w) e H^ (Q; D^) x H^{Q.;Du) is equal to the set of solutions of (3) . 

Proof. The identity 

{o{u)SI(p-y(p, w) = -{a{u](f>V(p,'Vw} + <cr(u)Vg0- V0, w) 

follows from Lemma 1 in [17]. The cut-off functional maybe used because of Lemma 1 above. The larger 
space of test functions does not change the set of weak solutions due to density and does not lead to 
infinite terms in (4). D 

We define the space X:= H^ (Q; D^) x H^ (fl; Du) and the afftne mapping 
and the nonlinear mapping 



<CToV^,Vl//> 

Vw) + {kv,w)r 



J) 



and the functional 



{[aiv)-ao)V(p,Vxfr) 



{a{v) \(p]V(p,Vw) - {a{v)Vgff>-V(p, w) 
\{h,w)R, 



where (p = g^ + cp and v = gu+v. Then equation (4) is in operator form 

Lx + Nx = b 

with X = {if), i2) e X. 

For {(f), i)), [y/, w) & X one has 

Na(p,v)AVf,w)) 

< CT°-(l|V(^||i2(n) + ||Vg^||i2(n)) 

• [II Vi^lli2(n) + (g° - go)ll Vw||i2(n) + Wyg^hnn) II i^WiHa)) 

< (l+||V(p||i2(n,R3))||{l^,M;)||x, 

where we use || (•, •) ||x to denote the natural norm in the product space X, in this case 



(5) 



(6) 



1(1//, w) I 



2 II ii2 
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Throughout the text we adopt the notational convention that for a function b one understands b = 
b - g(^ if b is a Greek letter and b = b - gt< if b is a Latin letter. We call H^ (O; D^) the first and H^{D.;Du) the 
second component of X. In this spirit we also refer, for example, to (cToV^, Vi//> as the first component of 
L. Furthermore, we distinguish between cf), which is a solution, and cp, which is a generic trial function. 



3 Existence and Convergence of Galerkin Approximations 

Consider a hierarchical family of subspaces {X„}„eM = {P„ x [/„}„efs^ whose union is dense in X. A Galerkin 
solution Xn e X„ of (5) is a solution of 

{Lxn + Nxn,y) = {b,y), VyeX„. (7) 

Lemma 3 examines continuity properties of L and N. 

Lemma 3. Let {yn}n = {i<Pn> Vn)}n be a sequence in X and y = {(p,v) e X such that (pn ^ (p,Vn ^ i> ^^ 
n ^ oo. Then Lyn -^ Ly weakly and Ny„ -^ Ny strongly in X* . 

Proof. Suppose there is a subsequence {v„(ic)}k and an e > such that 

llo'(z^(j(fc))V^„(fc)-o-(z^)Vg)||i2(n,K3)>£ VkeN. (8) 

The compactness of the embedding H^{D.) '-^ L^(0) and a corollary of the Riesz-Fischer theorem [21, 
p. 161] imply that there is a subsequence, also denoted {i',j(A:)}fc' which converges pointwise almost ev- 
erywhere. By possibly passing to another subsequence of indices we may also assume that {y(pn{k)}k 
converges pointwise almost everywhere. The sequence 

is bounded in each component by (o'°)^|V^„(fc)|^. From the dominated convergence theorem, in the 
form of (Royden, p. 270), it follows that the sequence {cr[Vn{k))'^(Pn(k)}k converges strongly in L^(Q,IR^). 
Observe that almost everywhere the poinwise limit of {(J{Vn[k))'^'Pn{k)}k is o'(f)V^, contradicting (8). 
Therefore Ly„ + Ny„ converges, indeed strongly, in the first component. It also follows that the terms 

criv„)\(p„]VcpneL^in,U^), 
a{v„)Vgff,-V(p„eL^{n), 
KVn^ L^idD.) 

converge strongly as n — ' oo. Hence {Ly„}„ converges weakly and {Nyn}n strongly to Ly and Ny in X* , 
respectively: 

,. {Ny„-Ny,Z) w ^* 1- /r r v ^ 

limsup = 0, VzeX : lim<Ly„ - Ly, z> = 0, 

" z^O \\z\\x " 

completing the proof. D 

The following lemma establishes a property ofL + N which is a variation of condition {S)o; a concept 
introduced by Browder, see [5] or [28, IIB, p.583] . 

Lemma 4. Let y„ = {{(pn, Vn)}n be a sequence in X and y=[(p,iJ)^X such that 

yn — y, (9) 

Lyn^Nyn-b, (10) 

lim<Iy„ + Ary„,((^„,0)> = <fo,((p,0)>, (11) 



lim<Ly,j + Wy„, (0, Vn)) = {b, (0, £>)>. (12) 

n 



Then yn^ y strongly. 



Proof. Adapting the argument of the proof of the previous lemma it follows analogously that {a ( y „) Vi//} n 
converges strongly in L^(r2,[R^). Using the strong convergence in (*), one obtains 

< limsup„<cr(f„)V((^-^„),V{^-(^„)> 

= limsup„[<a(z;„)V^ ,V^> -2<o-{y„)V(^„,V(^> + <CT(f„)V(^„,V^„>) 

= limsup„(<cr(f„)V^„,V^> -2(o-{y„)V(^„,V(^> + (o-(i^„)V(^„,V^„>) 

= limsup„{-{aiv„)y(pn,y(p) + {b,[(p,0)) + {aiv„)Vgfi,,Vi(pn-(p))] 

'=' -{b,{cp,0)} + {b,{cf>,0)} = 0. 

Therefore cpn converges strongly and Lemma 3 becomes available. Hence {Nyn, yn) -^ (Ny, y) and 

{b,y} = limn{Lyn + Nyn,y} =lim„<Ly„,y } + {Ny,y}, 
{b,y) = limn{Lyn + Nyn,y„} =lim„<Ly„,y„> + <Ary,y>. 

The weak continuity of L implies that lim„<Ly,2, y> = (Ly, y), cf. [10, p. 422]. Therefore 
= lim(Iy„, y,j - y> = lim<Ly,j, y> - 2<Iy,j, y> + {Lyn, yn) 

n n 

= lim<Ly, y> - 2<Ly„, y> + {Lyn, yn) = lim<i{y -yn),y- yn) • 

n n 

It follows from the coercivity of the linear part of L that y„ ^ y in X. D 

Let Tn'.Xn^ Xn be defined by y„ = ry„, where y„ = {(pn, Vn) ^ ^n is given as the solution to 

{Lyn + N[(pn,Vn),iW,w)) = {b,{M/,w)), {_i^,w)eX„ (13) 

with yn = (.(pn> Vn) ■ Algorithmically an iteration with T corresponds to a method with the primary variable 
Vfi and the dummy variable (p „ as (pn does not explicitly appear in the next step of the iteration. 

Lemma 5. There exists a radius r, independent ofn, such that the range ofTn belongs to 

Br:={yeX:\\y\\x<r} 
for aline N. 

Proof. Let y„ = [ipn, Vn) = Tny = Tn[0n, Vn)- The first component of (13) gives, with y/ = (pn, the identity 
= {a{Vn)y{(Pn + S<p)'^Vn)- Thus, with the above Poincare-Friedrichs inequality for 7/^^^; £></)): 

\Wn\\]ji(^) < {(y{Vn)^(Pn,^(f>n) = - {o {v n)^ g<p ,^ (p n) ■ 

The Cauchy-Schwarz inequality now gives \\(pn\\H^{p.) ^ \\^g4>\\L'^(Q.,'si^)- Recall (6) with (i//, w) = [(pn, Vn) 
and [(p, V) = {(pn, Vn)- The Cauchy-Schwarz inequality, (1) and the coercivity of the linear part of L give 
theboundedness of i>„. D 

Observe that the fixed points of Tn are exactly the Galerkin solutions in the sense of (7). 

Lemma 6. The mapping r„ has at least one fixed point x„ . 



Proof. We have that Tn : Br n Xn ^ B^ n Xn- The First Lemma of Strang (Braess, 2007, p. 106) implies 
that the Galerkin solution of a linear elliptic equation changes continuously in the H^ -norm as the diffu- 
sion coefficient is varied in the L°°-norm. Therefore, looking at the first component in (13), cpn depends 
continuously on Vn, taking the equivalence of norms in the finite -dimensional Xn into account. With 
(pn determined, Vn can be computed from {(p„, v„) = L~^{b- N{(pn, Vn))- Lemma 3 showed a sequen- 
tial continuity of property of L and iV which guarantees that the finite- dimensional Galerkin restrictions 
Xn -^ X,i axe continuous. Equally the Galerkin restriction Xn -^ Xn of the affine mapping L~^ is contin- 
uous. This means that Tn is a continuous map !„ : Bj- n Xn ^ Br n Xn, so Brouwer's fixed point theorem 
gives the existence of a fixed point x„. D 

It is a direct consequence that Galerkin solutions exist for all n e N and (that at least one of them) 
are contained in B, . The next theorem conceptually builds upon Proposition 27.4 in [28, vol. II B] where 
Lemma 4 is replaced by (S)o. 

Theorem 1. There exists a subsequence of Galerkin solutions {Xn(k)\k = ii^Mk), Vn(k))}k (^nd an x= {cp, v) 
in X such that Xn(k) ~^ ^ strongly in X and x solves (5). If the solution x of (5) is unique then the whole 
sequence converges. 

Proof It follows from (7) that <(I -i- N)Xn, y> -^ (b, y) for all fixed y e Xm, m e N. Observe that L^N i?, 
a bounded operator, see (6) for N. Thus with {x„}„ also the sequence {{L + N)Xn}n is bounded. Conse- 
quently, (L + N)Xn ^ fo in X* as n — ' oo, see Proposition 21.26 (c),(f) in [28, vol. II A] . The sequence {x„}„ 
is bounded in the reflexive Banach space X. Thus there exists an x e X and a subsequence {Xn{k)}k with 
•^n(fc) ^ X as A: ^ cxD. It follows from (7) that 

lim<Lx„(fc) + Arx„(fc),((^„(fc),0)> =lim<b,(^„(fc),0)> = <b,(^,0)>, 

k k 

lim<Lx„(fc) + Nxn(k]< (0, v„^k))} = lim<fo, (0, v„(k))} = (b, (0, v)). 

k k 

Therefore, it follows from Lemma 4 that Xn(k) -^ x as k ^ oo. The continuity established in Lemma 3 
yields (L + N)Xn(k) ^ {L + N)x as k ^ oo and therefore (L + N)x = b. If the exact solution is unique, then 
the Galerkin approximations can only have one accumulation point. D 

We remark that this convergence result applies, in particular, to conforming finite element methods. 

4 Regularity 

In this section we investigate how the regularity estimates for the Poisson problem with mixed boundary 
conditions, derived in [24] , carry over to the Joule heating problem. These bounds are sharp in the Pois- 
son setting. In general they are also sharp for the Joule problem — noting that equation (2) takes the form 
of Poisson's equation when choosing a constant a. 

The underlying question is whether additional regularity can be gained if the type of the boundary 
condtions only changes at re-entrant corners: 

(C) D is a creased domain. 

For the full definition of creased domains we refer to [24] ; here we only highlight that in the setting of 
the Joule heating problem the key conditions are that Q is a bounded Lipschitz domain, that D^ and Du 



are open and non-empty and that dD^p and dDu are not re-entrant, meaning that the angles between D^ 
and N(f, as well as between D„ and Ru are strictly less than Ji. 
Let "Ke c |r2 be the interior of the polygon with the vertices 

, (1,1), (l-e,l), |0,^ + e 



(0,0), (e,0), |l.--£ 



In the statement of the following theorem we let D e {Dff,,Du} and R e {Nff„Ru}. Also 1/p + l/p' = 1 
and -Bf'^ denotes the scale of Besov spaces. 

Lemma? ([24]). There exists an e = EidQ.,D,R) in (0, 1/2) such that Poisson's equation is well-posed in the 
spaces 



V eW'' , (Q), 



\v eiw!"' ,[n;D) 



2-s- 



P'P, 



v\d =geB^'^iD), dyV\R =heB';'_\{R], 



(14) 



whenever {s, lip) ^'K^. 

Consequently we assume for the remainder of this section that g^ £ Bf^{D^), gu e B^'^iDu) and 
h e B^i^^iRu). We also choose 



£:=mm{E[dn,D^,N^),e{dn,Du,Ru)}^iO,ll2). 
We aim to prove existence of solutions with the smoothness and Lebesgue indices 



2-£ 



s = 



t = 



8 + £ 

IT' 



p. 



l-£ 



q- 



12 

4-e' 



Then (s, l/p), {t,llq) e "Ke and 



W'' _, (Q) = W^-' (O) c C°'t (O) n Wf (Q) c C(n) n Wf (Q) 



(15) 



and 



W^ . (Q;D) = H^'«(Q;D)cL2^(n)cI^(n) 



(16) 



due to the Sobolev embedding theorem. Note that W^f ^ (f^) Ian = ^f ''' (3^2) has the dual space 

2 — 5 £ I c. 



cf. [24]. Also 



B^J^I^{dm = BP;_\{d^), 



W'^ , (Q) c W,2-' (Q) c Vl^."-^ (Q) = W'' , (Q). 



Indeed the embedding of W'^ j (O) into M^'' ^ [D.) is compact for £ e (0, 1/2) since 



s + - 



t + - 



l\ 1 



^/ 



-£. 



We set 



Y = W^^ 1 (Q; D<^) X VV^^^ ^ (Q; D J , Z = W^^ ^ (Q; D^) x W^^ ^ (O; D J . 



(17) 



Corollary 1. Let s, p and k e L"°{R) be as above. There exists an e = e{dQ.,D,R) in (0, 1/2) such that 
Poisson's equation is well-posed in the spaces 



(18) 



v\d eB^'^iD], KV + dyV\R eB^f^{R) 

whenever {s, lip) e JCg. 

Proof. A standard energy argument ensures that, given Dirichlet boundary conditions on D with data in 
B^'^iD) and Robin boundary conditions with data in B^^^(i?), the mixed Poisson problem has a unique 
solution V in H^^D.). This function v is also the unique solution of the system with Neumann boundary 
conditions 

h-KveLHR)<^ [Bfo^' iR)) * = [B^'Xo^R)] * = B^s-i ^^^ on R, 

giving the required regularity by Lemma 7. The result now follows since \\kv\\ gPf (^j ^ II ^ II h^ (n;D) • n 

Now, motivated by (2), we consider the following modified weak formulation for (0, u)eY 

{Vu,Vw) + {ku,w)r^^ = (a{u)V(f)-'SI<p,w) + {h,w)R^ Vic, J 

with i^f/, w) e Z. Notice that again (j) = (p-g^ and u=u-gu. 

Lemma 8. The set of solutions of (3) which belong to Y is equal to the set of solutions o/ (19) . Moreover, if 
{(p, u)eX solves (3) then it solves (19) for all [xf/, w] e W^in;D^) x W^in;Du). 

Proof. One only needs to consider the first equation of (3) and (19). Let {(f>, U) eV solve (19) . Let {de}e be 
an approximate identity and a£:=a{u) *de and ii/£:=y/ae. Then 

/ V0-Vi//dx'='lim / {——-■V(t)]if/dx, 
Jn t^oJnV ai; > 

as o E ^ o{u) in L°°(D). Now choosing !//■£ as test function 

lim / Vrf) • VWf dx = lim / f V(/)|i//-£dx=/ fVcr(u)- V(/)]wdx, 

but also 

lim / V0 • Vi//£ dx = lim / cr^V^- Vi//dx+ / fV</)- Vergil// dx 

= / o-(w)V(/)-Vi//dx+ / (Vo-(w)-V(/))v/dx. 
Jn Jn 

Subtraction shows that {(f>,ii) e Y solves (3). The other direction follows from re-arranging the above 
identities; the test spaces Z and W^ {Q; D^p) x W^ (Q; Da) have to be adapted to the size of trial function 
spaces. D 



It is convenient to define the operators 



whiere 



Si: Y 


-^ 


% . 


S2 : Z^ 


-^ 


% - 


S3: Y 


-^ 


w* , 


S4 : Zu 


— ' 


K - 


^f=K^- 


in;D^], 



{(p, v) '-' [w-' {a{v)\/ip-\/cp, w) + {h, w}rJ, 



-wP 



w^ = w:P__^{n;D^), 



2-s- 



Z,, = W;,{Q;DJ, 



p 



The notation indicates that operators map into dual spaces with the associated test functions if/ and v. 
Let / be the identity map and 



S:=(/,S4loS3)o(S2^oSi,/). 



(20) 



Given an initial pair {(p, v), Si returns in the first component the right-hand side of (19) and thus S2 ^ o Si 
gives an update of the first component. This, together with an unchanged v, is passed into S3 and then 
S^^ to return first an updated right-hand side and then an updated second component. 

Lemma 9. The operators S maps continuously into Z. 

Proof. Lemma 7 ensures that S2 ^ and SJ ^ map continuously into W'^ ^ [O.]. Notice that, due to (16), 



-——■V(p, a{v)\V(p\^ el— (Q) = i" (O) c wf ^ (Q) . 
a{v) V ; \^ 2-S-- j 

As heB^:![iRu],KeU^iRu), v e C(Q.) also Kv,he[wf_ _A^])*. Therefore Si and S3 map into W* and 



W* 



n 



Va 



Lemma 9 gives access to Schauder's fixed point argument provided -^ is not too large in relation to 
other parameters of the problem. Let Ci to C4 be the embedding constants of 



^2-s-^^^'^"^ 



L—e[n], 



wf.in] -^w^^.iD.], 



^+? 






wl^_Am\ 



tP 
'2-s- 



wr ,{nr, 



respectively By abuse of notation we denote by || Sg ^ II and || S4 ^ || the operator norms of the linear parts 
of S2 ^ and S4 ^; that is of S2 ^ and Sj ^ if g^ and g„ were 0. 



Lemma 10. There exists a positive constant 

C* =C*({C/}f^j,||S2^||,||S4 II, 0-°, II g0 II 5P.P, II giJgP.P, II h II ^P.p) 

such that whenever \\ ^ \\l°°(u) ^ C* then there is a ball B c F such that S maps B into B. 



(21) 
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Proof. By Holder's inequality 



2-s-p 

^Cil|2^||ioo(K) ||Vz;|| 12 llVmll 12 



Similarly, 



\\a{v)y<p-ycp\\ . = sup / ^ ^' (22) 

2-5- — 

We need to bound the first component /oS2^oS^ = 82^081 and thesecond component S^^oS3o(S2^oS^,/) 
of S, cf. (20). Suppose that {(p, v) are contained in the ball B = {y e Y : WyWy < r}. As S2 is an invertible 
affine function, there are generic constants C with a parameter dependence as indicated in (21) such 
that 

US^'oSyKip, v)\\^, („.. ) < C||2^lk«(R) r^ + C. (23) 

t+— 
1 

For the second component notice that the right-hand side of (22) only depends on <f> and not v 



\[Sl'oS3o{S2'°Si,im(p,v)\\^^, (^.^_^) <C(C II 2^ II i<^(K)r2 + C)Vc. 



owing to (23) 

Therefore || (S2 ^ o Si) [<p, v) \\ is bounded by a quartic polynomial of the form 

CM a' ||2 ,-4 I /^ II cr' m „2 , /-■ 

ll-5r||^„(K)r +C||-||L<»(K)r +C. 

At radii where it intersects the first diagonal r ^^ r the operator S maps B into B. The existence of such 
an intersection point is guaranteed if || ^ || l°°(u) is sufficiently small. D 

The above lemma is consistent with the analysis of the linear problem with a constant a in the sense 
that a suitable C* can always be found as 2- becomes small. We now turn to Schauder's fixed point 
theorem to obtain the existence of solutions in Z. 

Theorem 2. There exists a positive constant 

C* = C*({Ca-=pllS2 ^11, 1154^11, o-°,||g0|l5P,P(£,^,,||g„|l5P.P(£,_^),||/z|l5P.P(^j) 

such that whenever \\ ^ \\l,°°(K) ^ C* then (19) has a solution in Z. 

Proof. The compactness of S follows from Lemma 9 and (17). Now the result is a consequence of Lemma 
10 and Schauder's fixed point theorem. D 
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While the global regularity estimate in Y is sharp in the setting of creased domains (up to the distance 
oi{t,\l q) to the boundary of 'K^], as a comparison with the Poisson problem shows, more regularity is 
seen away from the boundary. For the next theorem we assume {ip, U) e Y; however, it is not relevant 
whether this is established with the above fixed-point argument or otherwise. 

Theorem 3. Let Qq be a relatively compact Lipschitz domain in Q.- Qq <s ^- Let [cp, U) eY bea solution of 
(19) . Then (p,ue W/ (Qq) for all s e (1, oo) . //cr e C~ (K) then (f),ueC"° [Th) . 

Proof. Let {OjI/gn and Qoo be smooth domains with Q, ^ 0.i+\ and D.i c Q^o <£ ^ for all i e N. Without 
loss of generality we may assume that the boundary data g^ and gu have extensions from the boundary 
onto Q such that g(p,gu £ C°°(Qoo)- Fix i eN\ {0}. Let d be a smooth function Qj -^ [0, 1] such that 
(ildfij =OandC;la,_i = 1- Then <;';0 and <;";« solve Poisson's problem on n, with homogeneous Dirichlet 
boundary conditions and the right-hand sides 

(i[l^yu■V(f>]-(i^g^-2yi;i■y(f>-4>M:i, 

(i[amV(p-y(f>]-(iAgu-2Vi;i-Vv-vA(i 
in Lt^iD.i), respectively. According to Theorem 9.15 in [14] a solution to Poisson's problem with ho- 

6 g 

mogeneous Dirichlet conditions belongs to Wg^ "(O/) if the right-hand side belongs to L4^(Qj) and 

6 12 

^ e (l,oo). Therefore (p,Ue W^^'^ [D.i) c W^"^'' [D.i). Substituting £ by 2e one may pass from i to J - 1 
and repeat the argument. We conclude via induction that (f>,u e W^{D.i) for any 5 e (l,oo), recalling 
that a negative ^^ corresponds to a right-hand side in L"°[D.i). Indeed 0, w e W^&i) for all i e N and 
5£ (l,oo). 

Now let a e C°°(IR). Theorem 9.19 in [14] states that if the right-hand side in Poisson's problem is in 
W^(n/) then the solution belongs to Wl^^&i) with k e {1,2,3,...} and 5 g (1,oo). Leibniz' rule if ■ g)^^^ = 

ly^o ( J )/'■'' •g**'"^' shows that for a given 5 e (1, oo) a s' e (1, oo) can be chosen such that ^^V0, o-(u)V0- 

V0 e W^{Q.i) if 0, w e W^^^{Q.i). Hence induction over k, coupled with a shift from Q,- to Q;-i as above 
to impose smooth boundary conditions, shows that (p,Ue W?(Qo) for all A; e I\l and s e (l,oo). Use of the 
Sobolev embedding theorem concludes the proof. D 

Remark 2. Assumption (C) is made to establish sufficient regularity of elliptic equations with Lemma 7 in 
the context of non-smooth domains and mixed boundary conditions. Also in other settings corresponding 
elliptic regularity results are available and the above analysis can be transferred with minor modifica- 
tions. We point for the pure Dirichlet problem to [19, 12], for the Neumann problem to [18, 12]. A related 
approach for the mixed problem on smooth, non-creased domains is proposed in [25]. 

5 A Priori and A Posteriori Error Analysis 

In this section we present a priori and a posteriori error bounds for conforming finite element approxi- 
mations. We first present in an abstract form that the Galerkin method is quasi-optimal and that the error 
X - Xn can be bounded using the dual norm of the residual of the approximation. In the second part we 
choose a particular approximation technique, namely a conforming /z-adaptive finite element method. 
We use interpolation estimates to bound the errors in terms of mesh size and polynomial degree. 
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5. 1 Abstract error bounds 

In the theorem we assume that (0, u) e Y. There is also an assumption on small data which relaxes as the 
coupling of the equations, measured by the Lipschitz constant of a denoted Cy , weakens. In order to get 
the correct dependency of Cy we introduce a scaling factor (t) of the second component of equation (5) 
in the proof of the theorem. Beside the Lipschitz constant Cy we let Cg be the embedding constant from 
H^ into L^, Cg be a Poincare-Friedrichs constant and 

C5(0) = C7C8(l + C9)||V0|li3(a)max(l,g°-go + ||Vg<^|li3(n)), (24) 

Ceicp) = ff°((l + C9)(C8l|V(/)||iS(n) + Cn\\yg<p\\LHn)^ + g°-go), 

We note that C5 is directly proportional to C7. 

Theorem 4. Suppose there is a solution x={(j},U)eY of equation (7) which satisfies 

(l-<5)2ffo 



C5 ((/)): 



C6((/)) + {l-5)CTo 



foraSe (0, 1) with the constants fi'om equation (24). Then the solution is unique. Furthermore, ifxn e Xn 
is its Galerkin approximation then the fijllowing a priori and a posteriori error bounds hold: 

\\x-x„\\x< inf \\x-yn\\x, 

\\X-Xn\\x < \\LXn + Nx„ - b\\x'. 

Proof. We introduce for {ip, v) eX the norm: 

We pick a solution x = {(p, U) e Y. For any function pair (!//„, Vn) = (i/^n + gcp, iJn + gu) with (i//„, Vn) e Xn, 
we have, 

\\[(p-(pn,U-Un)\\\T <{o{u„)y[(p-(pn),y{(p-(pn)) 



+ T^ll V(W - Un) ||^2jj^j + T^ll sficiu - Un) 11^2 



(Ru 



<<0-(W„)V((/) -(/)„), V((/) -!//«)> +T^<V(M-U„),V(W-y„)> 

+ T'^{K{U-Un),U-Vn)R,, + {Cf(_UnWi(p-(pn),yiV^n-(pn)) 
+ T^<V{M - U„),yiV„ - Un)) + T^{KiU - U„), Vn - U„}r^^ 

(*) a° 

< \\(_(p-(pn,U-Un) II xA\i(p-'^n,U- Vn) II X^ 

+ C7C8l|V(/)||i3(j^)||w-M,j||^i(n)||V((/)„-1//„)||i2(n) 

+ T2||o-(u)[01V0-O-(M„)[0„lV0„||i2(n)l|V(z;„-U„)||i2(n) 

+ CgT'^\\yg^\\r^(n)\\0{u)y(p-0{Un)y(pn\\mm\\^n-Un\\HHD.) 

a° 

< U(p-(p„,U-Un) II x^U4>-Vn,U- Vn) II X^ 

+ C6T2||V((/)-(/)„)||i2(n)||V(i;„-u„)||i2(n) 

+ C5(||V((/)„-1//„)||i2(n)+T2||V(y„-U„)||i2(n))l|V(M-U„) 11^2(0), 
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where we in (*) use that {a{u)\I(p,\I\f/) = {a{Un)y(pny^'H^) = for any ■^/ e X„. We now use the triangle 
inequality on II V(i^„ -(/)„) ||i2(n) < l|V(0- 0„)||i2(Q) + ||V(0-i//„)|li2(Q) and II V ( i;„ - uj 11^2(^5 < ||V(m- 
Un)\\L^(Q.) + \\^{-U-Vn)\\LHQ.) and introduce do = (T°a-^ + Cecr'^r + €5(7-^7-^ + CsT"^ to get, 



|(0-0,;,U-M,;)|||i <CioU<p-<Pn,U-Un)\\x^\\i(p-'Vn,U-Vn)\\x^ 
+ C6T^||V((/)-(/)„)||i2(o)||V(M-U„)||i2(Q) 

+ C5(||V((/)-(/)„)||i2(o) +T^||V(M- U„)||i2(o))||V(U- M„)||i2(n), 
Cio II {(p-(pn,U- Un) II x^ II {(f>-yfn,U- V„) \\ x^ 



< 



+ — (C5 + C6t2)||V((/)-(/)„)||2,j^j 



n'^^L'^iO.)- 



for any e > 0. We now want to find the maximum value of C5 while fulfilling ^(Cs + CeT^) < (1 - d)ao 
and [|(C5 + CeT^) + CsT^) < (1 -(5)t^, for some < 5 < 1, since we then can subtract (1 -6)\\{(p-(pn> u- 
Un)\\j^T on both sides of the equality sign. Algebraic manipulation reveals that by choosing e equal to 

e* = Ce+a-5)ao > ° ^"^ ^^ = ^^*^^ '= "''^q''"" > ^ l^ads to a value C5 = ^JriSk- U"^^*" *^ assumptions 
in the statement of the theorem we therefore have, 



lief) -cl)n,u-u„)\\ XT' <5 ^Cio(t*,o-o,o-'',C5,C6)||((/)-i//«,m- 



'nlWX'' 



for that fixed t* > 0. The a pr/orz bound now follows using the Poincare-Friedrichs inequality and simple 
algebraic manipulation, since the X^-norm and the X-norm are equivalent for fixed t. 

Suppose we have two solutions Xi and X2. The a pnonbound gives that any Galerkin approximation 
will eventually get arbitrary close to both solutions which means that they coincide, i.e. Xi = X2. 

Only the a posterioribound now remains to prove. It follows continuing ideas as in the a priori bound. 
We have, following the step (*) and below, 

||{(/)-(/)„,M-M„)|||-. <<0-(W„)V((/) -(/)„), V(0 -(/)„)> 

+ T^{y[U-U„),yiU-Un)}+T^{K[U-U„),iU-Un])R,, 
<{b- Ucpn, Un) - N{(pn, U„), {(p -(/)„, T^(U - Un))) 

+ C6T2||V{(/)-(/)„)||i2(n)||V(w-M„)lli2(n) 

+ C5l|V{M-U„)||i2(n)(l|V((/)-(/)„)||i2(n)+T^I|V(U-M„)||i2(n)), 

Under the assumptions of the theorem we can now repeat the exact same argument as above picking e* 
and T* in the same way. We get, 

II •^ - ^n III < II ^ - -^^n 11^.. ^ 5"^ max(l, (t* )^) {b - Lxn - Nx„,x - x„} 

<\\Lx„ + Nxn-b\\x'\\x-Xn\\x, 

completing the proof of the theorem. D 

Remarks. Since Cs{(p) andC%{(p) depend on l|V(/)||i3(n) abound on this normneeds to be computed in or- 
der to obtain a computable a posteriori error bound. In the setting of creased domains this can in principle 
be derived from the radius r in Lemma 10. 



14 



5.2 Error bounds for h -adaptive finite element approximations 

The local regularity result presented in Theorem 3 indicates that the problem is well suited for adaptive 
finite element methods. Suppose that T^ and T^ are decompositions of Q. into tetrahedrons. We let h t 
denote the diameter of element T and assume that T,^ and 7u are nondegenerate, i.e. there is a mesh- 
independent constant 7 such that 

max — <y, 
T dj 

where dj is the diameter of the largest ball contained in T. We denote the set of elements T' neighbour- 
ing r by Sr = {u T : r' n r 7^ 0}, assuming that elements are closed. Further let £/_^ , £/_^ be the set of all 
interior facets of the two meshes and let En^, £r^ be the set of boundary facets on A/^ and Ru, respec- 
tively. It is implicitly supposed here that A/^ and Ru are unions of edges. We denote the set of all elements 
neighboring e by Se = {ueiTnej^ 0}. The diameter of a facet e is denoted hg. The spaces of polynomials 
of total or partial degree less than or equal to k, defined on T, is denoted Pk(.T). Let 

P„:= {(peCin):(p\TeP,,{T]and(p\D^=0}, 1 
Un- {V eC[n]:v\TePg[T) and v\d„=0}, J 

and Xn = PnxUn- We let tt^ : H^{D.) -^ P„ and Uu ■ H^iO.) -^ Un be Scott-Zhang interpolants, as defined 
in [26] . We recall the interpolation bounds 

\\v-nuV\\wl(T) <^T~' l^lvv;(Sr). 1 

II II ^ui/2 I I r ^26) 

for v£. W2{0.;Du), re [l,fc-i-l], se [0,r] and where \-\w.' is the semi-norm only including the derivatives 
of order r. An analogous bound also holds for n^. From here on we drop the subscripts of ?!„ and n^p 
since it will be clear from the context which operator is meant. 

Theorem 5. Let x = {(f),u) e Y be a solution of (7) and x„ e X„ be its Galerkin approximation. Then we 
have the a priori error estimate, 

||r_r ||2 < V h2('-"Uft|2 + V h2(^-l)|/7|2 

with 1 < r < A; + 1 and \<s<£ + \, and the a posteriori error estimate, 

IU-^nll|< E h|T?r+ E ^er]\i+ E ^de.N 

TfiT^ eeE,^ beEn^ 

+ E h|Pr+ E ^epli+ E ^bpIr, 

reT„ ee£7„ e££fl„ 



where 



riTi(pn>Un) = 

r]e,li(pn<Un) = 

r]e,N{(pn, Un) = 

PTi(pn,U„) = 

Pe,d(pn,Un) = 
Pe,Ri(pn,Un) = 



y-a{Un)y(pn Ill2(r). 

V-[o{Un)y(pn\\\LHe)> 

V- a{UnW(pn WhHe)' 

y-iVUn + O-(U„)[0„lV(/)„) + O-(M„)Vg0-V(/)„||i2(y), 

V-([VM„]+Or(M,j)[0,jl[V(/)„])||i2(g), 

V-iVUn +a[U„)\4)„] V(/)„ )-HKM„-/l||i2(g). 
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Proof. From Theorem 4 we have, 

\\x-x„\\x< inf \\x-y„\\x^{\\ii>-n(t>\\lri + \\u-nufiY'^. 

The a priori pait of the theorem follows from (26) with s = 1. 

We turn to the a posteriorihound. For any y={cp,v)eX we have 

{Lx„ + Nx„-b,y} = {aiu„)y(l)n,y iq) - ncp)} + {V Un,y ii) - 71 D)} 

+ {KUn,V-nV) + {a[Un)\(f)n]y(l>n,yii>-7li»} 

-{a[Un)yg^-y(pn,v-nv)-{h,v-nv)R^^. 

We can subtract the interpolants because of Galerkin orthogonality. We apply Green's formula on the 
elements of the meshes together with the Cauchy-Schwarz inequality to obtain 

{Lx„ + Nx„-b,y)< Y, \\y-{a[Un)y4)„)\\L2^T)\\<P-^(PhHT) 

+ Y. \\v-l(r{Un)y(t>n]\\LHe)\\V>-^(PhHe) 

+ Y. W"^ ■i^Un + 0{Un)\(pnW<Pn))+Cfi-Un)yg^-y(pn\\LHT)\\^-T^^\\LHT) 
Te7u 

+ Y llv-([VM„]+cr(u„)[(^„l[V(/)„])||i2(g)||i)-7i:i;||i2(g) 

+ Y \\v-{yun + aiu„)\(f)„]y(l)„) + KU„-h\\iHg^\\v-ni)\\L2(^g^ 

Use of (26), with 5 = and r = 1, and the Cauchy-Schwarz inequality gives 

{LXn + Nx„-b,y) , ^ 2 2 V-u2 V-u2 U/2 

reT„ ee£7,„ ee£R„ 

The theorem follows by taking supremum over all yeX. D 

Remark 4. If the Dirichlet data g^p and gu are not traces of finite element functions then an extra data error 
term will appear in the error estimates. Similarly if exact quadrature is not used additional error terms 
have to be considered. We have neglected these standard terms in the analysis to make it more readable. 
For more details, see [6]. 

Remark 5. For additional flexibility one can consider hp finite element spaces. Theorem 3 indicates that 
even though the choice of the space Y is sharp, the interior regularity of the solution is higher This is 
an ideal setting for hp -finite element methods. The above a priori and a posteriori error bounds can be 
transferred to hp approximation spaces for which suitable interpolation operators are available. Some 
ideas for the construction of such operators are collected in [23] and references therein. 
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